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Abstract. Ice water path (IWP) and cloud top height (h t ) are 
two of the key variables in determining cloud radiative and 
thermodynamical properties in climate models. Large uncer- 
tainty remains among IWP measurements from satellite sen- 
sors, in large part due to the assumptions made for cloud 
microphysics in these retrievals. In this study, we develop a 
fast algorithm to retrieve IWP from the 157, 183.3 ± 3 and 
190.3 GHz radiances of the Microwave Humidity Sounder 
(MHS) such that the MHS cloud ice retrieval is consistent 
with CloudSat IWP measurements. This retrieval is obtained 
by constraining the empirical forward models between col- 
located and coincident measurements of CloudSat IWP and 
MHS cloud-induced radiance depression (T c i r ) at these chan- 
nels. The empirical forward model is represented by a look- 
up table (LUT) of 7' C i r -IWP relationships as a function of h t 
and the frequency channel. With h { simultaneously retrieved, 
the IWP is found to be more accurate. The useful range of 
the MHS IWP retrieval is between 0.5 and lOkgm -2 , and 
agrees well with CloudSat in terms of the normalized prob- 
ability density function (PDF). Compared to the empirical 
model, current operational radiative transfer models (RTMs) 
still have significant uncertainties in characterizing the ob- 
served 7’Cir-IWP relationships. Therefore, the empirical LUT 
method developed here remains an effective approach to re- 
trieving ice cloud properties from the MHS-like microwave 
channels. 


1 Introduction 

Ice clouds have profound impacts on the global energy bud- 
get (Stephens et ah, 1990), hydrological cycle (Chahine, 
1992), atmospheric structure (Ramaswamy and Ramanathan, 
1989) and circulation (Richter and Rasch, 2008). Cloud ice 
water amount is one of the largest sources of uncertainty 
in quantifying cloud-climate feedbacks and sensitivities. For 
example, the mean cloud ice water path (IWP) ranges from 
10 to 120 gm -2 in the tropics among a variety of global cli- 
mate models (GCMs) in the most recent 20th century Cou- 
pled Model Intercomparison Project Phase 5 (CMIP5) runs 
(Li et ah, 2012). Accurate cloud IWP measurements are criti- 
cally needed to guide model developments and reduce model 
uncertainties. 

However, observations of cloud ice have not met the 
requirement by climate models, showing several folds of 
IWP differences among various techniques (Wu et ah, 2009; 
Eliasson et ah, 2011). Until cross-instrument consistency is 
achieved, current cloud ice observations will allow too much 
variation in cloud properties and become insufficient for con- 
straining the model physics (Waliser et ah, 2009; Li et ah, 
2012). Difficulties for accurate IWP and microphysical mea- 
surements arise mainly from remote sensing in the presence 
of cloud inhomogeneity and sensitivity limitations associ- 
ated with each technique. On one hand, large spatial and 
temporal variabilities in cloud microphysics make it diffi- 
cult to compare ground-based measurements with remote 
sensing observations (Waliser et al., 2009). Hence, statisti- 
cal representations of cloud microphysics are assumed or pa- 
rameterized in order to enable satellite remote sensing (e.g., 
McFarquhar and Heymsfield, 1997). Even for simple opti- 
cally thin cloud, there are still a great number of uncertainties 
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in the assumptions made for the IWP retrieval. On the other 
hand, passive satellite sensors have limited penetration capa- 
bility to observe thick and dense ice clouds from space. As a 
result, only partial columns of IWP (pIWP) can be measured 
by passive sensors, and the column bottom varies with at- 
mospheric absorption, cloud amount, droplet size and phase, 
and cloud top height. These uncertainties about cloud col- 
umn create additional errors in the IWP retrieval using pas- 
sive sensors. 

As an active sensor, CloudSat radar provides an unprece- 
dented opportunity to measure the ice water content (IWC) 
profile and its vertical integral (i.e., IWP) globally since 
2006. The CloudSat cloud ice retrieval still depends on the 
cloud microphysics constrained by in situ and ground-based 
observations (Austin et ah, 2009). CloudSat data are confined 
in a narrow curtain (~ 1 km width) along the orbital track, 
and thus are used mostly for climatological and case stud- 
ies. Like other A-train sun-synchronous satellites, it samples 
only two local solar times (01:30 and 13:30LST) of the cloud 
diurnal cycle. However, CloudSat data still provide the best 
characterization of the vertical distribution of global cloud 
ice (Eliasson et ah, 2011), and can be used to cross-calibrate 
other techniques, especially the passive sensors with limited 
vertical resolution. (Wu et ah, 2009). 

Passive nadir-viewing microwave techniques such as Ad- 
vanced Microwave Sounding Unit-B (AMSU-B) and the Mi- 
crowave Humidity Sensor (MHS) have an advantage over in- 
frared/visible sensors in penetrating deeper into cloud layers 
to measure IWP. More importantly, MHS has a swath width 
of ~ 2300 km to capture synoptic-mesoscale systems in mo- 
tion as well as variabilities not captured from the curtain- 
only sampling by CloudSat. Instead of slicing a single ver- 
tical cross section of a hurricane, the entire cyclonic struc- 
ture can be mapped out with one MHS orbit. Since 1998, 
satellites carrying instruments like AMSU-B and MHS have 
been operational and now fly across the Equator at more than 
eight local solar times every day, the mosaic of which can be 
used for cloud diurnal cycle studies. Moreover, at microwave 
frequencies ice scattering signals are approximately linearly 
proportional to cloud ice amount in the path, resulting in 
a relatively straightforward relationship between IWP and 
cloud-induced radiance depression (Wu et ah, 2009). These 
advantages make nadir-viewing microwave sensors attractive 
for monitoring global long-term IWP. 

Retrieval of IWP requires radiative transfer models 
(RTMs) or forward models that relate cloud ice to the mea- 
sured radiance. The cloud ice models can be formulated ei- 
ther theoretically or empirically. RTMs are also widely used 
in climate models but primarily for calculating clear-sky 
radiative forcing from atmosphere gas, cloud, aerosol and 
surface. Although studies demonstrated the use of RTMs 
for IWP retrievals from AMSU-B/MHS channels, consid- 
erable uncertainties exist with RTMs in representing com- 
plex physical processes (e.g., land surface radiative fluxes, 
ice particle shape) and with oversimplified assumptions (e.g., 


plane-parellel atmosphere and cloud layers, cloud droplet 
size distribution, etc.). The errors in liquid drop size, sur- 
face emission/scattering, cloud layer height, and water va- 
por amount can all degrade the quality of the retrieved IWP. 
For example, the current operational IWP retrieval algorithm 
from Microwave Surface and Precipitation Products System 
(MSPPS), which is based upon a two-stream approximated 
radiative model solution (Zhao and Weng, 2002) at AMSU- 
B 89 and 150 GHz channels, was found to under-estimate 
IWP in comparison with other observations (Wu et ah, 2009; 
Waliser et ah, 2009; Eliasson et ah, 2011; Chen et ah, 2011). 
Contamination of cloud ice retrievals was also found over 
snowy/icy surfaces (Wu et ah, 2009). 

While further improvements are still needed for ice scat- 
tering calculation in the microwave RTMs, empirical forward 
models have been used for cloud retrievals (Holl et ah, 2010). 
Empirical approaches establish some ad hoc relationships 
between cloud ice variables and radiance/reflectivity mea- 
surements from the data themselves. Such empirical forward 
models are developed from a finite ensemble of observations, 
and are therefore limited to specific conditions, environments 
and dynamic ranges of the cloud variable of interest. The al- 
gorithms are usually fast in the form of a look-up table (LUT) 
and bypass the complex microphysical calculation in cloudy- 
sky radiative transfer in individual cases. Empirical methods 
have also been used in surface remote sensing where land 
properties are too complicated to be modeled or validated 
(e.g., Pulliainen and Hallikainen, 2001). 

In this paper, we develop an empirical model and retrieval 
algorithms for IWP using cloud-induced radiance depression 
(T clr ) from MHS at 157, 183. 3±3 and 190.3 GHz. The empir- 
ical forward model is obtained by regressing MHS Tcir radi- 
ances on collocated CloudSat IWP and cloud top height mea- 
surements in the tropics. The sequential estimation method is 
then used to retrieve IWP for all MHS footprints. The instru- 
ments and methodology will be described in Sect. 2, followed 
by the detailed retrieval algorithm in Sect. 3. An evaluation 
of the retrieved products and the associated errors is going to 
be given in Sect. 4, followed by the summary in Sect. 5. 


2 Instruments, data, and methods 

2.1 Description of data sets and models 

The data sets used in this study are Level- 1 brightness tem- 
perature (7 b) from MHS, ice water content (IWC) from 
CloudSat, and Modern Era Retrospective-Analysis for Re- 
search and Applications (MERRA) three-hourly analysis 
variables on a 1.25° x 1.25° latitude-longitude grid. The two 
radiative transfer models used in this study are Joint Center 
for Satellite Data Assimilation (JCSDA) Community RTM 
(CRTM) and an ice scattering cloud radiance model (CRM). 
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2.1.1 MHS Tb, IWP and historical issues 

MHS is a cross-track scanning radiometer aboard the Na- 
tional Oceanic and Atmospheric Administration (NOAA) 
satellite 18, 19, European Organisation for the Exploitation 
of Meteorological Satellites (EUMETSAT) Metop-A, and 
Metop-B, which is a slightly improved version of AMSU- 
B onboard NOAA- 15, 16, and 17. MHS makes 90 footprints 
(with a beam width of 1.1°) continuously in each cross-track 
scan and the outermost scan angle is ±48.95° from nadir. 
For NOAA- 18, the MHS scan and satellite orbital altitude 
produce a nadir footprint size of 16 km at half-power field of 
view (FOV) and a swath width of 2200 km. The FOV size and 
swath vary slightly among satellites due to different orbital 
altitudes. MHS has five microwave channels, which are 89, 
157, 183.3 ±1, 183.3 ±3 and 190.3 GHz (for AMSU-B, the 
second and last channels are 150 and 183.3 ± 7 GHz, respec- 
tively). For consistency, these channels are labeled as CH#1- 
CH#5 hereafter. MHS CH#1, CH#2 and CH#5 are vertically 
polarized, and the other two are horizontally polarized (for 
AMSU-B, all five channels are vertically polarized). The 
designed radiometric noises (NEAT) for CH#1-CH#5 are 
0.22, 0.34, 0.51, 0.40, and 0.46 K, respectively (John et ah, 
2012). The 89 and 157 GHz are window channels, and those 
around 183.3 GHz water vapor absorption line are designed 
to profile the atmospheric water vapor. Under clear-sky con- 
ditions, the peak sensitivity of these 183.3 GHz channels oc- 
curs in the upper, middle and lower tropospheres, respec- 
tively. NOAA-15, 16, 17, 18 and 19 orbits drift slowly with 
time, while Metop-A and Metop-B are maintained at a sun- 
synchronous orbit with fixed Equator passing time (EPT). 

For ice particle scattering measurement, the higher- 
frequency channels (157, 183.3 and 190.3 GHz) work bet- 
ter for IWP retrievals because the Mie scattering is propor- 
tional to frequency to the fourth power. Scattering-based mi- 
crowave cloud remote sensing has some unique properties as 
well as limitations. First, it penetrates deeper into ice clouds 
than IR and visible techniques for cloud ice measurements, 
but can become saturated for very optically thick clouds (Seo 
and Liu, 2006; Arriaga, 2000). In the case of saturation, only 
partial cloud ice column pIWP can be retrieved. As shown 
in Seo and Liu (2006), the window channels near 183.3 GHz 
can penetrate a cloud layer with IWP as large as lOkgm -2 , 
which covers most of the IWP values observed by Cloud- 
Sat. However, in the case of graupels, or frontal astrostratus 
clouds, saturation may occur (Arriaga, 2000). Saturation is 
also more prominent in the oblique views than nadir (where 
the line-of-sight path is longer). 

Secondly, among all MHS/AMSU-B channels, CH#3 is 
most sensitive to water vapor because it is adjacent to 
the 183.3 GHz water vapor absorption line. The absorption 
from upper tropospheric water vapor, so-called “water va- 
por screening”, prevents CH#3 from seeing the surface and 
clouds in the lower troposphere. To some extent, CH#4 has 
quite amount of water vapor screening and can observe some 


ice clouds, but remains little contaminated by the surface un- 
less at dry, high latitudes. In other words, CH#4 can be used 
to distinguish between surface and clouds in the situation 
where other channels have difficulties, as will be shown in 
Sect. 2.1.3. 

Lastly, microwave radiances are dependent on scan angle 
at these frequencies. Under clear-sky conditions, the radi- 
ance may decrease with scan angle from nadir, as a function 
of the cosine of angle, due to the increasing path length at 
line of sight (LOS). This is similar to the 6.7 pm IR chan- 
nel where the longer LOS path gives a weighting function at 
a higher altitude, or cold temperature (Soden, 1998). Under 
the cloudy-sky condition, the radiance scan dependence may 
vary with cloud inhomogeneity as cloud size and distribution 
are often not homogeneous. In addition to the atmosphere- 
induced scan angle dependence, there are some instrument 
errors in all five channels that are scan dependent and asym- 
metric about nadir. These instrumental errors can severely 
degrade the quality of the retrieved IWP if not properly cor- 
rected. For example, there was a radio-frequency interference 
(RFI) problem in CH#3 and CH#4 of AMSU-B (Atkinson, 
2001; Buehler et al., 2005), and gain variations/degredations 
are found in CH#3-CH#5 of AMSU-B on NOAA-16 and 
NOAA-17 (John et al., 2013). MHS exhibits smaller scan- 
dependent biases than AMSU-B, but suspicious behaviors 
have been reported for CH#3 on NOAA- 18 and NOAA- 19 
and Metop-A (John et al., 2013). The MHS instruments on 
NOAA- 18 and Metop-A have so far shown the best overall 
radiometric calibration for all five channels. Since NOAA- 
18 has the closest EPT with CloudSat, it is used in this study 
to develop the cloud ice retrieval constrained by CloudSat. 
The radiances from CH#3 are not used because they are rel- 
atively noisier and provide little information on cloud ice. 
As in the main weather prediction centers, we use the Ad- 
vanced Television and Infrared Observational Satellite Op- 
erational Vertical Sounder (ATOVS) and the Advanced Very 
High Resolution Radiometer ( AVHRR) pre-processing pack- 
age (AAPP, v7) developed by Numerical Weather Prediction 
Satellite Application Facilities (NWP SAF) to process the 
LIB radiance data to obtain the further quality-controlled 
and calibrated L1C data. In the NOAA-18 MHS L1C data 
we have not found any systematic instrumental error. 

Weng et al. (2003) developed an algorithm to retrieve the 
IWP using ice scattering at 89 and 150 GHz, which is known 
as the NOAA operational IWP product. Their retrieval algo- 
rithm yields effective ice particle size and IWP with cloud top 
and base temperatures derived from simultaneous AMSU- 
A channels. A considerable fraction of false cloud detection 
was found with this method, mostly over icy/snowy surfaces 
and on elevated topography (Wu et al., 2009). The NOAA 
IWP has been reported to have significantly low values com- 
pared with radar and IR measurements (Holl et al., 2010; 
Eliasson et al., 2011). As an extended product, rain rate is 
derived from the retrieved IWP with an empirical polynomial 
relationship (Ferraro, 2007). The operational NOAA IWP 
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data, now integrated into the MSPPS in the CLASS website, 
will also be used in this study for comparisons. 

2.1.2 CloudSat IWC 

Launched into the A-train in April 2006, CloudSat has a 
94 GHz cloud profiling radar (CPR) to provide continuous 
cloud profiles along its nadir track. A CPR FOV size is 

1.3 x 1.7km. The cloud ice water content (IWC) product 
from 2B-CWC-RO (R04) is used in this study, which as- 
sumes a gamma size distribution of cloud ice particles. The 
CloudSat IWC retrieval is limited when the temperature is 
above 0°C; so is the liquid water content (LWC) retrieval 
at temperatures below —20 °C. Between 0 and —20 °C, IWC 
and LWC are retrieved separately and linearly interpolated to 
the intermediate temperature range (details of the algorithm 
can be found in Austin et al., 2009). Thus, large uncertainties 
are expected for this mixed-phase cloud regime, and/or in 
the ice cloud cases with large snow/graupel particles present. 
The vertical resolution of the IWC profile is 250 m. In our 
study, we interpolate it vertically to an evenly spaced grid 
(250 m resolution), and integrate the IWC between surface 
and 1 9 km to compute the total IWP. We also integrate the 
IWC profile from different bottom heights to represent the pl- 
WPs measured by MHS channels better. Compared with Holl 
et al. (2010), who used the CloudSat total column IWP prod- 
uct, our IWC integration approach is more meaningful for 
comparison with pIWP seen from MHS water vapor chan- 
nels, although the pIWP value is calculated on a profile-by- 
profile basis. Hereafter, we use IWP as the abbreviation of 
pIWP in our study to represent the MHS cloud ice column. 

CloudSat IWC has been validated with in situ, ground- 
based and other satellite IWC measurements (e.g., Austin 
et al., 2009; Wu et al., 2009; Protat et al., 2009). The un- 
certainty is claimed to be up to 40 % (Austin et al., 2009), 
which is much smaller than the divergences among various 
satellites and models, the latter of which often exceed 100% 
(Waliser et al., 2009; Eliasson et al., 2011). In this study, we 
treat CloudSat IWP as the “truth” to constrain the retrieved 
MHS IWP difference relative to that of CloudSat. Moreover, 
since microwave penetrates much deeper into ice clouds than 
IR/V1S channels, we expect our CloudSat-constrained algo- 
rithm to yield a better retrieval at large IWP values. 

2.1.3 Radiative transfer models (RTMs) and 
computation of 7” c i r 

The first step in cloud ice retrieval is to determine ice cloud- 
induced brightness temperature 7' C]r from raw radiance mea- 
surements (Wu et al., 2009, 2014). In this study, 7’ ar (cloud- 
induced radiance) is defined as the difference between the 
measured radiance, Tg, and modeled clear-sky background 
(also called “cloud-cleared radiance”), T ccr : 

Lcir = 7g — T CCI . (1) 
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Tdr also serves as a critical variable for cloud detection 
since every measurement has an uncertainty that may lead 
to a false alarm. T c ; r error is largely affected by uncertainty 
in the estimated Leer- Various methods have been developed 
to improve the accuracy of Leer estimation. Generally speak- 
ing, T ccr can be obtained using statistical differences between 
cloudy and clear skies (Wu et al., 2005), or using the radia- 
tive transfer model to estimate the clear-sky background from 
the current atmospheric state. Here we use the second ap- 
proach with the best estimate of local atmospheric state vari- 
ables (e.g., temperature, pressure, water vapor, ozone) and 
surface conditions (e.g., surface temperature, surface type) 
from a MERRA three-hourly assimilation data set from in- 
terpolation of adjacent grid points and closest local time. We 
allow relative humidity to exceed 1 00 % in computing clear- 
sky radiation. We also used the MERRA six-hourly finer-grid 
analysis product and European Re -Analysis Interim (ERA- 
Interim) data, but no statistically significant difference in 7' ccr 
distribution is found among the results so far in the tropics 
and subtropics. 

The JCSDA CRTM v2.0.5 model is employed to calcu- 
late T cct . CRTM is a fast radiative transfer model that uses 
an advanced doubling-adding method (Liu and Weng, 2006) 
to compute the radiances and radiance Jacobians at the top 
of the atmosphere for various instruments, with wavelengths 
ranging from visible to submilimeter. It includes scattering 
calculations for cloud, aerosol, gas molecules and surface if 
specified. As the key backbone of data assimilation (DA) sys- 
tems, the CRTM has incorporated most space-borne instru- 
ment information (e.g., spectral frequency, filter shape, and 
scan pattern), including AMSU-B and MHS. Therefore, it is 
also our objective to calibrate our cloud ice retrieval with this 
widely used CRTM for T cct estimation so that the IWP out- 
puts can be ready for the DA applications. 

Figure 1 presents the probability density functions (PDFs) 
of 7’ c , r , T CCI and 7fi from a month’s worth of MHS nadir mea- 
surements in the tropics. Warmer 7g values are mostly from 
the clear sky or surface, while colder 7g are the cases of ice 
clouds or snow/icy surfaces at a high elevation. The 7fi PDFs 
all have a broad peak with a standard deviation (er) that is 
so wide that the empirical 3cr cloud detection method (i.e., 
L’Bpcak - < 0 for cloud detection) used by many previous 

studies does not work well when applied directly to the 7g 
data (e.g., McNally et al., 2006; Gong and Wu, 2013). On 
the other hand, the 7’ c , r PDFs have a smaller standard de- 
viation because the CRTM-derived T CCI from MERRA data 
has removed a lot of clear-sky variability (Fig. la). The long 
PDF tail in the negative Lcir values is a distribution of cloudy 
radiances. Ideally, a perfect Leer model with perfect clear- 
sky input would produce a singular peak in the 7’ cir PDF at 
0 K, and all negative values that are smaller than the radiance 
noise would be classified as clouds. The uncertainty of T ar 
measurements, close to a Gaussian distribution, is reflected 
in the PDF spread near zero, especially in the positive half 
of the PDF. The ability to separate between cloudy and clear 
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Figure 1 . Probability density functions of T c j r (a), T ccr (b) and 7b (c) for CH#2 (black), CH#3 (blue), CH#4 (green) and CH#5 (red). Samples 
are from the NOAA-18 MHS nadir view during August 2010 in the tropics ([25° S, 25° N]). 


radiances is characterized by this standard deviation er, which 
can be computed from this portion of PDF. 

However, the CRTM does not always improve cloud de- 
tection. For example, the width of the CH#2 T ar PDF is not 
much narrower than that of 7 b, indicating limited skills of 
the CRTM in capturing the clear-sky variability. Large error 
in the calculated Tccr is found over mountains and arid areas, 
where it remains challenging for the CRTM to model sur- 
face contributions at CH#2. When excluding all land cases, 
the CH#2 T'cir can produce a PDF with a narrower width 
around zero (not shown). On average, the r c ; r error is ~ 5 K 
(one standard deviation), although it may vary from 7.5 to 
10 K. In the cloud ice retrieval later on, the generic value of 
5 K is used for all channels. In addition to r c ; r standard de- 
viation, we also calculate T c „ bias for each MHS channel 
since the clear-sky PDF should peak at zero. We find that the 
CRTM has a cold bias (~ 2K) at 157 GHz (see Appendix B 
for details), whereas the bias is negligible for other channels 
(Fig. 6). 

Moreover, Fig. la also reveals the dynamic range and pen- 
etration depth of the four MHS channels in measuring cloud 
ice. Ch#2 penetrates deepest into clouds. Benefitting from its 
low frequency (157 GHz) at which cloud scattering and wa- 
ter vapor absorption is lowest among the MHS channels, it 
produces the longest cloud PDF tail (black line in Fig. la). 
On the other end, CH#3 has the most absorption from wa- 
ter vapor, showing the smallest T’cir dynamic range. It has a 
slightly broader distribution in the positive half of T’cir PDF, 
compared to those of CH#4 and CH#5, indicating that ei- 
ther the upper-tropospheric water vapor from MERRA or the 
CRTM calculation at CH#3 contains greater uncertainty. 

Two popular operational RTMs are also used to explore 
how the observed T c i r -IWP relationship is simulated by 
models. These two RTMs are CRTM, and a multi-stream 
“cloudy-sky radiance model” (CRM) that is currently used 
by the Microwave Limb Sounder (MLS) team to retrieve ice 
cloud properties (Wu and Jiang, 2004). They both show a 
certain lack of ability to capture the observed T c i r -IWP re- 
lationship, which results in large biases in data assimilation 
and increases the uncertainty of IWP retrieval (Appendix C). 


2.2 Collocated and coincident MHS-CloudSat 
measurements 

Collocated and coincident measurements (collocations here- 
after for briefness) are the incidences where two or more 
sensors observe the same location at the same time. These 
measurements provide useful pairs for instrument calibration 
(e.g., John et al., 2012), cross-validation of a particular vari- 
able (e.g., Wang et ah, 2010), or development of new retrieval 
methods (e.g., Lamquin et ah, 2008). In this paper, we will 
be focusing on the last application. 

The requirements for collocated-coincident measure- 
ments may vary, depending on the variability of the specific 
variable. Since most of the atmospheric state variables (e.g., 
wind, temperature, humidity) change relatively slowly and 
continuously with space and time compared to fast processes 
like clouds, their requirements for collocation and coinci- 
dence should be a bit more relaxed and the allowed windows 
for space and time should be consistent. In other words, the 
uncertainty of collocation due to spatial variations should be 
comparable to one of coincidence due to temporal variations. 
Another factor in defining the requirements for collocation 
and coincidence is to assure enough samples for statistics. 
For the A-train sensors, sample size is usually not a problem. 
On the other hand, such a near-perfect collocation is rare 
between radiosonde and Global Positioning System (GPS) 
measurements (Sun et ah, 2010). Neither occurs frequently 
for two satellites that run in different orbits. Adjustment of 
the collocating criteria becomes necessary and important in 
these situations. 

In this study we use NOAA-18 measurements to find 
collocated-coincident cases with CloudSat because NOAA- 
18 has the closest LST to CloudSat orbit among all opera- 
tional satellites with the MHS/AMSU-B instruments (Holl 
et ah, 2010). The requirements for collocation and coinci- 
dence are 10 km in space and 15 min in time, which yield a 
total of 6 x 10 5 MHS samples in the tropics (25° S and 25° N) 
from June 2006 to March 2011. Holl et ah (2010) obtained 
one order of magnitude more collocated NOAA18-CloudSat 
measurements with a requirement of the same time differ- 
ence but 15 km in distance, of which the number increase is 
roughly proportional to the area differences between the two 
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Figure 2. Total number of collocated and coincident MHS foot- 
prints as a function of scan angle between June 2006 and 
March 2011. 

distance criteria. The sensitivity of the retrieval algorithm to 
the choice of collocation criteria will be discussed in the next 
section. 

Because of the close orbits between NOAA-18 and Cloud- 
Sat, the number of collocated measurements peaks at the 
MHS nadir angle and drops off similarly at the left and 
right view angles (Fig. 2). There is no significant scan angle- 
dependent sampling bias, which would be a factor to consider 
in the derived Tcir —IWP relationship. The number of collo- 
cations decreases sharply at oblique views with scan angle 
6 > 35°, which may affect the statistical significance of the 
derived Tcir— IWP relationship. 

In the case of highly inhomogeneous clouds, larger uncer- 
tainty is expected for the IWP within MHS FOV, as CloudSat 
footprints cover at most 6.7 % of the area of an MHS foot- 
print. As a matter of fact, multiple CloudSat cloud profiles 
often correspond to an MHS footprint because the CloudSat 
footprint (~ 1.5 km) is much smaller than the spatial range 
of the defined collocation. Thus, we average all the CloudSat 
IWP values within the collocated MHS FOV to represent the 
mean IWP for the MHS footprint. The same procedure is ap- 
plied to calculate the mean cloud top height (/z t ) at that MHS 
footprint, where each individual CloudSat A t is obtained by 
searching for the highest level where IWC > lOmgm -3 . 

3 Empirical T c -, r model and IWP retrieval 

3.1 Empirical r c j r -IWP relationships 

For nadir-viewing sensors like MHS/AMSU-B, negative T c i r 
is caused primarily by ice cloud scattering instead ofby emis- 
sion. Mie theory shows that T ar is proportional to cloud IWP 
at microwave wavelengths and to the fourth power of fre- 
quency. As the ice cloud becomes radiatively thick, cloud 
self-extinction prevents T c [ r from penetrating deeper to sense 
the entire IWP, but rather, it has a sensitivity to pIWP Hence, 
an empirical T c i r — IWP relationship is derived in the follow- 
ing format: 

r cir = r cir0 .(i- e - iwp ^), (2) 


where 7' clr o is the coldest 7’ C]r (i.e., saturation value) and H 
is the parameter to determine where T c { r becomes saturated. 
Both r c j r o and H depend on frequency and can vary with 
cloud top height (/z t ), instrument view angle, and temperature 
lapse rate (y ) in the upper troposphere. In this study, since we 
focus on the tropical region where the lapse rate variation is 
small, these parameters are assumed to be only a function 
of channel frequency and cloud top height. For small IWP 
values, r c i r ~ 7’ clr o(—IWP//7 ), which is a linear relationship 
as described by Wu and Jiang (2004) for Aura MLS. As also 
suggested by Wu and Jiang (2004), H could be a function of 
cloud profile shape and the ice-to-water mixing ratio inside 
clouds, but these dependencies have secondary effects on the 
7’cir-IWP relationship. 

To derive the empirical 7’ c]r -IWP relationship, we first sort 
all collocated measurements, CloudSat IWP (averaged onto 
MHS footprints) and MHS T c ; r at near-nadir views (scan an- 
gle e [—5°, 5°]), to generate a joint PDF separately for each 
MHS channel. As shown in Fig. 3, the 7’ c j l -IWP relation- 
ships are scattered with the PDF peaks in good agreement 
with Eq. (2). We then fit the 2-D PDF to obtain 7’ cn -o and H 
parameters in Eq. (2), which is the solid curve in Fig. 3. The 
fitting is carried out as follows: (1) to determine T^ro from 
the coldest T c „. We search all 2010 MHS nadir data and the 
coldest 7^^ as 7’ C |,-o for each channel. (2) We then compute 
H in Eq. (2) with the ordinary least squares method by fit- 
ting the 7' clr and IWP values at peak 2-D PDF (black dots in 
Fig. 3) using the derived from step 1. 

The fitted curves represent bulk characteristics of the joint 
PDF. Compared to a linear fit, the residual variances de- 
creased by at least 50%. However, the joint PDF of CH#3, 
showing a steeper relationship for ?cir and IWP at colder fir 
values, is not represented well by Eq. (2). Moreover, r c ; r PDF 
becomes flat at small IWP values (IWP < 0.5 kgm -2 ), indi- 
cating the lower limit of T c „ sensitivity to IWP. The spread 
of 2-D PDF reflects both natural variability and collocation 
error of the 7' C]r -IWP relationship. One of the cloud vari- 
abilities that affect the 7’ CH -IWP relationship is the cloud top 
height (At). 

To examine the dependency of H on At, we further 
sort the collocated measurements into three height groups 
using the mean A t (defined as the highest altitude at 
which IWC reaches lOmgm -2 ) computed from CloudSat 
cloud profiles: 9.5 < At < 10.5 km, 1 1.5 < At < 12.5 km and 

13.5 < A t < 14.5 km, each group separated by 1 km to avoid 
overlapping of the regression lines. These three height groups 
account for about 48 % of all near-nadir collocated mea- 
surements. We then apply the same fitting procedure to ob- 
tain 7^0 an d H for each height group as in Eq. (2) (solid 
thick lines in Fig. 4). For the three height groups, a cloud 
bottom height (Ab) is calculated to be within 7.5 ± 1.6 km, 

8.5 ± 2.5 km and 7.4 ± 2.5 km, respectively, so that the col- 
lection represents tall and thick deep convective clouds in the 
tropics. The measurements with A t below 8.5 km and above 
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Figure 3. Contours of 2-D PDFs (normalized by the maximum 
value of the histogram) of collocated CloudSat IWP (abscissa, av- 
eraged onto MHS footprints) and MHS r c j r (ordinate) for CH#2- 
CH#5 at near-nadir views (scan angle e [—5°, 5°]). Black dots lo- 
cate the peak of the 2-D PDFs, and the thick solid lines are regres- 
sion curves using Eq. (2). 


15.5 km are too few to obtain a statistically robust T’cir-IWP 
relationship. 

Fig. 4 shows that T c „ is more sensitive to IWP for clouds 
with higher h t , except for CH#3, where the situation is re- 
versed. This variation in 7’ clr sensitivity is expected, accord- 
ing to the sensitivity expression from a conceptual cloud scat- 
tering model (Eq. 6.3 in Wu and Jiang, 2004) 


7 c i r 

^ Tscat J A I i , 

r ceff 


( 3 ) 


where r ce ff is the cloud effective optical depth that is posi- 
tively correlated with IWP, 7’ scat is the cloud scattering radi- 
ance from a convolution of the upwelling and downwelling 
radiation, and 7 ab is the background clear-sky radiance be- 
neath the cloud. For a given channel. Tab remains the same 
no matter how thick or thin the ice cloud is. For thick, high /i t 
clouds, 7’scat is colder due to more contributions from higher 
altitudes, resulting in a larger T cn sensitivity to IWP. The 
CRM used by Wu and Jiang (2004) for Aura MLS predicts a 
similar but weaker h t dependence, due to the fact that MLS is 


Table 1. Look-up table for the parameters of the joint I WP—h\ re- 
trieval. 



^cirO 

[K] 

H) 

[kgm- 2 ] 

Cl 

[gm -3 ] 

c 2 

[mgm -4 ] 

CH#2 

-172 

21.45 

-1.9875 

0.5625 

CH#4 

-140 

17.021 

-0.4078 

0 

CH#5 

-155 

29.6511 

-2.26214 

0.038156 


a limb sounder that observes the cloud side in its LOS rather 
than the cloud top seen from a nadir sensor. 

Since we do not have an accurate model of H dependence 
on h t , a quadratic function is assumed to interpolate and ex- 
trapolate H(li t ) to the cases beyond the values at the ob- 
served ht, i.e., /it = 10, 12 and 14 km. The coefficients in 
Eq. (4) are solved from the observed H values for three /i t 
groups with the mean values at 10, 12 and 14 km. Including 
7’cirO, all the parameters of the empirically derived r c i r — IWP 
relationships for CH#2, CH#4 and CH#5 are listed in the 
look-up table 1 . 


H = co + ci/? t + C2/i t 2 (4) 

The /q-dependent H parameter allows a simultaneous re- 
trieval of /q and IWP. By including or constraining /i t in the 
retrieval, it improves the IWP retrieval. Other approaches, 
e.g., using IR channels from the CO 2 slicing method (Kahn 
et ah, 2008) may be used in the future to constrain /j t in the 
IWP retrieval. As seen in Fig. 4, the error bar for each of 
the three cloud groups is smaller than one without the height 
separation. Relaxing the collocation requirements would in- 
crease the number of measurements for statistics, but we find 
that it does not reduce the error bar of the derived 7 c i r — IWP 
relationship. 

To complete the empirical model for the 7 c j r — IWP rela- 
tionship, we need to extend the parameters listed in Table 1 
from the near-nadir case to all MHS scan angles. For off- 
nadir views, to account for longer off-nadir LOS (f is the 
local zenith angle), the ^cir o£f _ nadir nee ds to be multiplied by 
cost; to achieve an equivalent nadir Tcirr assuming plane- 
parallel cloud layers. This is not a bad assumption in the case 
where clouds are not opaque. For opaque clouds, inhomo- 
geneity plays a more important role in relating off-nadir and 
nadir views. In other words, the scan-angle correction for 7’ cir 
is a function of 7’ cir as well. Thus, we develop an empirical 
solution for this correction, which is given in Appendix A. 

3.2 Joint retrieval of IWP and /i t 

The IWP and /i t are retrieved using the sequential estima- 
tion approach as described in Rodgers (2000) and Livesey et 
al. (2006). Equation (5) in Livesey et al. (2006) is quoted as 
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Figure 4. PDF peaks (uncertainties given as error bars) and the cor- 
responding regression lines based on Eqs. (2) and (4) for clouds 
with h\ between 13.5 and 14.5 km (black), 11.5 and 12.5 km (blue) 
and 9.5 and 10.5 km (red) for CH#2-CH#5 at near-nadir views. 


Eq. (5) below: 

x ( 9 +i) = x (q) + Sx [K T S;\y - T^)) + S ~\a - x ( ^)] (5) 

(q) annotates the r/th step of iteration. In our case, x = 
[IWP, h t \ is the retrieved result, y = [T c ;„ 2 , 7’ CII - 4 , 7’ c j r , ] is the 
observation, and 1 can be calculated using Eq. (2) and 
x^\ K is the Jacobian matrix, which is defined as 


r 3 T’cir i 


r Cm e -iw?/H -i 

H ' e 

SIWP 


97cir 

L dll t J 


WIWP. e -IWP ,H. (ci+ 2 c 2 h t) 


Plotted in Fig. 5 are the analytical solutions of K using 
the coefficients listed in Table 1 and Eq. (6). 97’ c j r /9IWP 
(left column of Fig. 5) monotonically increases with IWP 
for all three channels without any singularity point or mul- 
tiple solutions. Flowever, 37' clr /9/z t (right column of Fig. 5) 
has a singularity point at h t = 18 km for CH#2, where mul- 
tiple solutions exist. For CH#5, multiple solutions can also 
occur for /z t . if we define the bottom of the <)T cu -/dhi curve 
as /q critical , then the smaller the IWP is, the higher /7 tcntlca , is. 
For instance, /q critical at IWP = 3.0 kg m -2 is 18 km, meaning 
that if the cloud has li t > 18 km and IWP = 3.0 kg m -2 , the 









Figure 5. Analytical solutions of the two components of the Jaco- 
bian matrix K: 3T’ c j r /9IWP with fixed lit (left) and dT^/dht with 
fixed IWP (right). For the left column, the fixed h\ value increases 
from 6 (thin, blue) to 20 km (thick, red) with an interval of 2 km. 
For the right column, the fixed IWP value increases from 0.5 (thin, 
blue) to 1 8 kg m -2 (thick, red) with an interval of 2.5 kg m -2 . 


retrieved /i t has a possibility of being underestimated. The K 
matrix responses at CH#2 and CH#5 suggest that the /z t re- 
trieval could significantly underestimate the truth when cloud 
top is above 18 km, especially for thick, dense clouds. 

a — [IWPo, h, Q ] is the a priori (initial guess) of x. In prac- 
tice, if T’cir from all three channels is less than —5 K, there is a 
strong possibility of ice cloud presence, and h to is set to 5 km 
to speed up the convergence of the iteration. Otherwise, h to 
is set to 0 km instead. The initial guess of IWPo is always set 
to 0. Once the iteration begins, a is forced to equal to x iq ] to 
avoid “artificial preference” of retrievals to the a priori. That 
is to say, the last term on the right-hand side of Eq. (5) can be 
eliminated. The total number of iteration steps is set to 20 re- 
gardless of whether the final results converge or not. Within 
each iteration, IWP <9) is not allowed to exceed 25kgm -2 or 
become negative, and the h[ q) value must be within the range 
of [0, 18] km. The lower bounds assure physically meaning- 
ful solutions. The upper bound of /z t is where CH#2 and 
CH#5 are problematic in retrieving a trustable h t with the set 
of coefficients listed in Table 1. Therefore, the protection of 
the ht solution again significantly under-evaluates the li t for 
those high, dense clouds. Nevertheless, IWP rarely exceeds 
25 kgm -2 , and the monochromaticity of K with respect to 
IWP assures the robustness of IWP retrievals. 

S y , S a and S x are the matrices describing the error covari- 
ances associated with the measurements, the a priori, and the 
final retrieval results, respectively. S y = [5,5, 5] 2 K 2 as the 
measurement error is estimated to be 5 K (Sect. 2.1.3). S a as- 
sociates with the CloudSat IWC retrieval error, which is es- 
timated to be less than ~ 50 % (Austin et al., 2009). In prac- 
tice however, S a defines the step allowed to jump in each 
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iteration, which needs to be small in very nonlinear cases 
where multiple solutions exist and large steps could result in 
an unstable retrieval. Since the retrieval function is mono- 
tonic for all channels, a large step S a = [6kgm -2 , 6 km] 2 is 
chosen, as in the so-called Newtonian iteration, to acceler- 
ate the retrieval convergence. Once S y and S a are fixed, S x at 
each iteration step can then be calculated from Eq. (5), which 
is shown as Eq. (7): 

S I = [K r S- 1 K + S- 1 r 1 . (7) 

The retrieval is not carried out if Tea at all three chan- 
nels is greater than 5 K, a strong indication of a clear sky. 
In that case, we directly assign a clear-sky flag to the scene. 
CH#2 radiance is excluded for retrievals over arid areas be- 
cause of its contamination by surface signals. This is real- 
ized by checking land pixels with Tea < — 5 K for all three 
channels (i.e., ice cloud likely). As long as this criterion is 
not satisfied, only CH#4 and CH#5 are used for the retrieval 
over land, whereas CH#2 is always used over oceans. Af- 
ter retrieval, the IWP value that has the standard deviation 
(V5JTT) greater than or equal to itself is flagged as “bad 
quality”; so is h t . The rest is flagged as “good quality”. 

4 Assessment of IWP and h t retrievals 

Comparisons of IWP retrievals have been challenging and 
sometimes even confusing because not all sensors measure 
the same portion of pIWP. Different cloud bottom and top 
heights can affect the cloud ice sensitivity and retrieval re- 
sults. For MHS, the channel penetration depth varies with 
water vapor loading above cloud and with liquid water 
amount inside clouds if it is a mixed-phase case. In addition, 
cloud inhomogeneity along LOS introduces more uncertain- 
ties to this comparison task. Active microwave sensors such 
as CloudSat do not have the penetration depth issue for most 
clouds. In this study we treat its IWP as the truth when com- 
paring it with the measurements from passive sensors (e.g., 
Wu et al., 2009). Since the retrieval algorithm developed 
here is constrained by CloudSat IWP, the IWP retrieved from 
MHS is expected to be statistically close to CloudSat cloud 
ice. In other words, MHS penetration depth is “extrapolated” 
to reveal the total column IWP using this algorithm. In this 
section we compare the PDFs of monthly IWP as well as the 
mean IWP maps for MHS and CloudSat data. 

4.1 Comparison of IWP PDFs 

Normalized PDF has been used to compare the cloud ice 
products and sensitivities from multiple sensors (Su et al., 
2009; Wu et al., 2009). The fundamental assumption of this 
approach is that cloud ice should have the same probability 
distribution if both sensors are measuring the same ensem- 
ble of clouds (e.g., in similar latitude regions and local time). 
Unlike the apple-to-apple comparison, vast data can be di- 
gested in one PDF plot that reveals ample information. The 
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Figure 6. PDFs of CloudSat IWP (grey thick line; smoothed over 
15 CloudSat footprints and integrated between 5 and 19 km), all re- 
trieved MHS IWP (black solid line; from all views), retrieved MHS 
IWP that is quality controlled (black dots), and MSPPS IWP (black 
crosses; from all views; from NOAA-18 only) for August 2010 in 
the tropics. 


basic philosophy of this approach is that the variable of in- 
terest should have the same probability of observing a certain 
value with what nature shows within the product’s visibility 
range. Therefore, if the probability is smaller (greater) than 
that from the truth, the variable (e.g., cloud occurring fre- 
quency) is under-(over-)estimated. The PDF comparison also 
overcomes the instrument geometry difference, as explained 
in Wu et al. (2009). 

As expected for the CloudSat-constrained retrieval, MHS 
IWP PDF agrees well with CloudSat, as shown by the grey 
and black lines in Fig. 6. The decreasing probability with 
IWP reflects the natural variability of cloud ice. CloudSat 
IWPs here are 15-FOV averaged values in order to mimic 
the MHS footprint diameter, which is slightly steeper than 
the original (non-averaged) PDF, or a higher (lower) possi- 
bility at smaller (larger) IWP. The averaging effect (< 10 % 
in PDF values) is negligible compared to the differences 
among various data sets/retrievals. When all good and bad 
retrievals from the 90 MHS views are included, the PDF 
(solid black line) in Fig. 6 rises more sharply at small IWPs 
(~ 500 gm -2 ) due to the arbitrary retrieval suppression for 
negative IWP values and false detection of clear-sky scenes. 
The dropping PDF at IWP < 500 g m -2 is mostly noise. 
When the quality flag is applied to exclude bad retrievals, the 
PDF (dots) agrees better with CloudSat at IWP > 300 gm -2 . 
At large values (> 8 x 10 3 gm -2 ), our algorithm tends to 
slightly over-estimate IWP when compared to CloudSat. 

The PDF of NOAA operational MSPPS data (crosses in 
Fig. 6) is lower than CloudSat at all IWP values. At large 
IWP values (IWP > 10 3 gm -2 ), it differs by 10 times or 
more, indicating that the operational product significantly 
underestimates cloud ice, compared to CloudSat. This low 
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Figure 7. 7 cir at CH#2 (a), CH#4 (b) and CH#5 (c) and retrieved NOAA-18 MHS IWP (d), h t (e) and MSPPS IWP (f) for Hurricane Earl 
at 01:54 LST on 31 August 2010 (Cuba is the island to the left of the plot). The IWP (ht) calculated from the collocated and coincident 
CloudSat overpasses (averaged onto MHS footprints) are marked by color crosses that share the same color bars with MHS. Note that the 
blank in the IWP map means that the IWP retrieval is not performed or has failed because it is below the sensitivity level of this algorithm, 
while CloudSat overpasses show a zero IWP value in most places away from the hurricane. 


bias was also reported in other studies (e.g., Waliser et al., 
2009; Eliasson et al., 2011). 

The quality of our cloud ice retrieval is demonstrated 
clearly in a scene over Hurricane Earl on 31 August 2010 
(Fig. 7a-c). 7’ c i, in all three MHS channels captured the struc- 
ture of Hurricane Earl very well, showing the eye, eye wall 
and spiral rain bands. CH#2 radiances, penetrating the deep- 
est, reveal more ice cloud structures than other channels. The 
retrieved IWP from our algorithm (Fig. 7d) retains most of 
the fine structures in CH#2 and also shows a hint of an ad- 
ditional two outer arms. Although the values of these arms 
are below the noise level, they are probably real because 
they are also present in the geostationary satellite 1R im- 
age (not shown). The IWP from CloudSat overpasses (col- 
ored crosses) have slightly larger values than MHS, whereas 
MSPPS operational IWP (Fig. 7f) are significantly smaller 
than CloudSat and our retrievals. This hurricane case also 
highlights the value of MHS IWP in studying the 2-D at- 
mospheric dynamics and cloud structures that are not cap- 
tured by the CloudSat curtain sampling. Using the CloudSat- 
constrained IWP measurements, we can obtain good spatial 
and temporal coverage from the MHS/AMSU-B sensors on- 
board all operational satellites. 

Retrieved h t (Fig. 7e) also agree reasonably well with 
CloudSat, especially at the hurricane periphery and the eye 
wall, but are lower by ~ 4 km over the hurricane deck (13 km 
versus > 18 km). This is probably due to the fact that the 
cloud top at the deck is dominantly higher than that at the 
hurricane periphery, i.e., higher than 18 km. They hence ex- 
ceed the upper limit of the reliable h t retrieval range from 


our algorithm. Pixel-by-pixel comparisons are done for some 
other cases that have CloudSat cloud tops lower than 18 km, 
and the h t retrieval seems quite promising (not shown). Nev- 
ertheless, the ht retrieval here is mainly to improve IWP re- 
trieval, rather than the purpose of scientific study. 

As the first CloudSat-calibrated column-wise IWP mea- 
surement that has excellent spatial coverage, the MHS IWP 
retains numerous potential usages for model input, for val- 
idation of other instrument measurements and for model- 
observation comparisons in the future. 

4.2 Geographic distribution of IWP 

Monthly mean IWP maps show good correlation between 
MHS and CloudSat cloud ice for August 2010 (Fig. 8), where 
the correlation is 0.81 in the tropics. Sampling error is evi- 
dent in these maps. With a relatively coarse grid box (5° x 5°), 
the CloudSat monthly maps (Fig. 8a, c) are spotty due to a 
lack of swath coverage. This sampling is also aligned to the 
westward-traveling fast cloud systems, leading to cloud ice 
spikes (e.g., eastern Pacific) and scatters (e.g., Amazon rain- 
forest) on the CloudSat maps (e.g., Amiridis et al., 2013). 
The sampling bias is largely mitigated by the 90-FOV MHS 
swath, maps of which look much smoother instead (Fig. 8b, 
d) because it overpasses one grid box 6 times as often as 
CloudSat on average. If the footprint size is taken into con- 
sideration, MHS could pass every corner of each 5° x 5° grid 
box in the tropics by as many as 42 times within a month, 
while CloudSat covers only 4 % of the area in the tropics. The 
major features between the two data sets agree well, espe- 
cially in deep convective regions where IWPs are large. The 
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Figure 8. Monthly averaged IWP from CloudSat and MHS ascending (a, b) and descending (c, d) orbits during August 2010. MHS IWP is 
averaged over all views. Data are sampled to 5° x 5° grid boxes. 


day-night differences in ice cloud thickness seen in CloudSat 
are also evident in the MHS maps, e.g., in central America 
and central Africa. 

Interestingly, in the scatter plot of MHS and CloudSat 
IWPs at a logarithm scale, the correlation is not along the 
1 : 1 line, showing a higher bias in MHS at smaller IWP 
values. The overall regression yields IWPcioudSat = (0.83 ± 
0.017)IWPmhs ~ 14.7 [gm“ 2 ], shown as the blue dots in 
Fig. 9b. The — 14.7 gm -2 offset partly comes from elevated 
topographies, e.g., the Andes, and from deserts, e.g., central 
Australia. The bias is slightly worse during night (MHS de- 
scending orbit) than during the day (MHS ascending orbit). 
If CH#2 is included for the MHS IWP retrievals over land, 
the high bias would increase over Australia, which may sug- 
gest a warm bias in MERRA surface temperature or aprob- 
lem with surface emmisivity in that region during nighttime 
(i.e., a cold bias of T cir for CH#2). It is suggested that our re- 
trieval algorithm has some limitations over complicated sur- 
face conditions, which will be discussed in the next section. 
Part of the — 14.7 g m -2 offset is caused by the fact that MHS 
tends to slightly over-estimate IWP with respect to Cloud- 
Sat, especially for thick and dense clouds. Besides, CloudSat 
probably misses some convections due to its sampling bias, 
for instance, over the Amazon rainforest and the maritime 
continent. Visual comparison between MODIS ice cloud op- 
tical depth (Fig. 7 of Meyer et ah, 2007) and MHS IWP 
shows better agreement in these regions. It is worth mention- 
ing that collocated MHS-CloudSat retrieved IWP showed the 
same feature. However, by subtracting the square root of the 
error (S x ), MHS IWP does not have such a positive bias (not 
shown). Therefore, despite the fact that MHS is still noisy 
below 0.5 kgm -2 at a single retrieval, the error estimation is 
very reasonable, and helps in filtering out the bad retrievals. 

4.3 Limitations of the algorithm 

The empirical forward model (7’ tn -IWP relationships) and 
retrieval algorithm presented in this paper are designed for 
tropical regions and have difficulties in retrieving IWP over 
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Figure 9. (a) is the same as Fig. 6, except from different latitude 
bins (see sub-titles for the latitude range), (b) is the scatter plot 
of MHS (abscissa) and CloudSat (ordinate) gridded monthly mean 
IWP for latitude bins between [25° S, 25° N] (blue filled dots), [25°, 
30°] N, and S (light blue triangle). The map grid size is 5° x 5°, and 
data are then smoothed by a two-point window along latitude and 
longitude before making the scatter plot. The black thick (thin) line 
marks the 1 : 1 (1 : 5 and 5 : 1) ratio. 

elevated topography and desert. In the cases of mixed-phase 
clouds or excessive water vapor abundance above cloud tops, 
the retrieval error for IWP might increase. The major causes 
of the biases over land are likely the CRTM surface emissiv- 
ity error in modeling the surface radiation, or surface temper- 
ature error in the mountain and desert regions in the MERRA 
data. Since CH#2 radiance contains surface signals (CH#5 
sees arid and snow surfaces as well), uncertainties in surface 
temperature and emissivity will induce T c „ biases. As a mat- 
ter of fact, we do see a systematic warm bias of 2 K in CH#2 
T c ir (Appendix B), which could be due to the instrument cal- 
ibration error or T c cr model error. Moreover, the PDF of T CCI 
for CH#2 (Fig. lb) extends to a temperature as low as 220 K, 
which is strong evidence of contamination from clouds or 
cold surfaces (i.e., ice pack on mountains can also cause this 
low 7b). With the development of a neural network approach, 
the initial guess of T CCI could be used to improve the cloud 
ice retrieval over complicated surface conditions (Chen and 
Staelin, 2003). 
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(a) Retrieved IWP Error 




if MERRA water vapor is too dry or too wet above clouds. 
The water vapor impact was only evaluated using CRM with 
different water vapor profiles, assuming variability within the 
uncertainty of observed upper-troposphere water vapor. The 
water vapor impact is found to be small and negligible in 
these CRM simulations (less than 5 % with doubling water 
vapor amount above clouds). 


Figure 10. (a) PDF of relative error (VSiwp/IWP, %) of one- 
month retrieved IWP between 30° S and 30° N. Contours are in 
log scale. The peak of the PDF is roughly marked by hand to rep- 
resent an exponential decrease (thick dashed line), (b) Collocated 
CloudSat-MHS IWP retrievals for July 2010 within the same lati- 
tude band. The thick solid line is the 1 : 1 line, and the other two 
thin solid lines mark the 1 : 2 and 2 : 1 ratios, respectively. 


The parameters in Table 1 assume that the atmospheric 
temperature lapse rate y is constant in the tropics. As pre- 
dicted by Wu and Jiang (2004) using CRM, the 7’ C „-IWP 
relationship is also a function of y (Fig. 6.10 therein). Eval- 
uating the PDFs of retrieved IWP outside the tropics, we 
find that the PDF of extratropical IWP starts to oscillate at 
its large-value tail (Fig. 9a) compared with that in the trop- 
ics (Fig. 6). The quality-controlled PDF in this case is still 
comparable with CloudSat PDF though for this bin assum- 
ing a 100% tolerance level of the PDF difference. Hence, 
our algorithm is expected to perform well within latitudes of 
30° N, S, but degrades in the extratropics. The mean vertical 
temperature profiles are also similar to those in the tropics up 
to 30° in latitude (Fig. 6.9 in Wu and Jiang, 2004). At middle 
to high latitudes beyond 30°, the quality-controlled PDFs are 
too low or even alter its shape, and the retrieved MHS IWP 
merely correlate with CloudSat IWP (not shown). In future 
algorithm development, y should be treated as an indepen- 
dent variable, such that the algorithm can be applied for IWP 
retrievals at higher latitudes. 

Liquid clouds occur frequently below 5 km, where temper- 
ature is usually greater than 0°C (Riedi et al., 2001), which 
may have little impact on CH#4 but can significantly af- 
fect CH#2 and CH#5 7’ cir . For deep convective clouds, liquid 
droplets can be lifted to a much higher altitude. The mixing 
of liquid droplets into ice cloud enhances the cloud emission 
contribution at microwave frequencies and hence decreases 
the Tcir sensitivity to IWP. Wu and Jiang (2004) showed that 
this impact could be as large as 30-50 % in a strong mixed- 
phase case, which alters the relationship in Eq. (2) with dif- 
ferent parameters. Therefore, mixed-phase clouds can con- 
tribute significantly to the spread of the 2-D PDF shown in 
Figs. 3 and 4. 

Water vapor above and inside cloud plays a screening role 
in reducing the sensitivity to IWP, in a way similar to liquid 
droplets. Since CH#4 and CH#5 are water vapor channels, 
they are sensitive to the water vapor abundance above and 
inside ice clouds. As a result, T c cr calculation could be biased 


4.4 Error analysis 

The retrieval error (S x ) provides a direct estimate of the re- 
trieval uncertainty. This uncertainty is completely indepen- 
dent of CloudSat IWP retrieval uncertainty. There are fur- 
ther sources of error, for example, the imperfect “coincident 
collocation”, the uncertainty induced by limited collocation 
samples, cloud misclassifications, etc. The total error is a 
combination of the three. This section will be focused on de- 
lineating the error sources and quantifying them one by one. 

The retrieval error S x is a dependent of the observational 
uncertainties from MHS (,S' V ) and from the forward model 
(K matrix), represented by Eq. (7). For small IWP values 
(thin ice cloud), S x is dominated by S y (i.e., MHS instru- 
ment noise). Since S v is fixed at 5 K for all three channels, 
this algorithm loses its sensitivity for IWP values smaller 
than 0.5 kgm -2 , which is evident in Figs. 9a and in 10b. For 
large IWP values (thick or precipitating ice cloud), S x is con- 
trolled by the uncertainty of the forward model. In the case 
when 7),,- saturates (i.e., smaller than 7’ cll -o in Table 1), the 
forward model uncertainty is infinite as the Jacobian curve 
of the K matrix with respect to IWP flattens out from 7’ ar <) 
onward. However, the K matrix is still very sensitive to IWP 
as long as Tc ir does not saturate, as shown in Fig. 5. There- 
fore, the induced error in S x from the forward model uncer- 
tainty is relatively small for large IWP values. This is con- 
firmed by Fig. 10a, where the percentage range of retrieval 
error ( V Stwp /IWP) decreases quickly from above 100% at 
small IWP values to as low as 20% at large IWP values 
(IWP > 1 kgm 2 ). Error is also reflected in the retrieval re- 
sults, shown in one month’s worth of collocated CloudSat- 
MHS IWP retrievals in Fig. 10b. It is apparent that collocated 
retrievals are highly agreeable above 1 kgm -2 , while MHS 
tends to overestimate thin cloud IWP. 

Ice cloud misclassification is an unavoidable issue for any 
cloud retrieval technique. Cloud misclassification of this re- 
trieval algorithm is partly induced by the beam-filling ef- 
fect and mismatch of CloudSat and MHS footprints spa- 
tially and temporarily. It is hard to separate these two ef- 
fects, as they are essentially the same. They are the major 
cause of the spread of the 2-D PDF and the uncertainty bars 
in Fig. 3 and Fig. 4, respectively. As a CloudSat footprint 
is much smaller than an MHS footprint, cloud inhomogene- 
ity within an MHS FOV cannot be captured fully by averag- 
ing CloudSat footprints within the corresponding MHS FOV. 
Since we further relax the collocation and coincidence crite- 
ria, mismatches also occur. Both of them result in the cloud 
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misclassification as described in the above paragraph. How- 
ever, imperfect clear-sky background ( T CCI ) can also cause a 
bias in 7’ C]r , which can also induce cloud misclassification. 
One thing to notice is that the error of clear-sky radiance is 
internally included in Sy, as S y is directly estimated from 
T’B-T’ccr- 

As truth is not given in our case, CloudSat is again used 
as the “truth” since the philosophy of this paper is to align 
MHS measurements with those of CloudSat. Two types of 
cloud misclassification exist: for one thing, MHS treats a 
footprint as an ice cloud being present, but collocated Cloud- 
Sat measurements averaged onto the same MHS footprint do 
not report a positive IWP value. This is named “Type I” mis- 
classification; for the other, which is on the contrary, it is 
called “Type II” misclassification. All collocated CloudSat- 
MHS observations between 30° S and 30° N are compared 
through the entire year of 2010 to generate the statistics. For 
Type I misclassification, there is an 18.1 % chance that MHS 
will detect an ice cloud but CloudSat will not, among which 
only 1 % of the MHS retrieval reports an IWP greater than 
0.5 kg m -2 (noise level). That means that using our tech- 
nique, MHS nearly never misclassifies a clear-sky scene as 
cloudy sky and reports an IWP value beyond its noise level. 
Interestingly, the 18.1% misclassification cases do not pre- 
fer mountain or snowy regions like the Andes, albeit there is 
a slight enhancement of occurring frequency over the Aus- 
tralian desert. The latter is expected as discussed in Sect. 4.3. 
For Type II misclassification, there is a 45.8% chance that 
CloudSat will detect an ice cloud but MHS will not. Among 
these cases, 90.5 % of the time CloudSat reports an IWP 
value less than 0.5 kgm -2 . All in all, our technique shows 
a very strong confidence level for retrieved values greater 
than the detection threshold (2 and 4 % misclassification rate 
for Type 1 and II, respectively), which means that our tech- 
nique nearly never misses a thick ice cloud, even when the 
surface signal is complicated. On the other hand, our tech- 
nique misses or gives large retrieval uncertainties of thin ice 
clouds like cirrus, which are not very important hydrologi- 
cally (Fig. 8) but are important to the radiation budget. 

It is difficult to justify the magnitude of other errors, 
which are blended in the error bars of the fitted curves in 
Fig. 4. The IWP error bar can be as large as ±1 kgm -2 at 
IWP = 3 kg m -2 for CH#4. Given the fact that the retrieval 
error should be smaller than any of the error bars generated 
from each individual channel, it is fairly reasonable to claim 
that the total retrieval error would be smaller than 100% 
when the retrieved IWP is beyond the detection threshold. 

The overall retrieval error from this algorithm is quite 
small. As a matter of fact, it reflects the “precision” of the 
instrument and the forward model rather than the total er- 
ror. One should always keep in mind that this is the error 
on top of the CloudSat IWP retrieval uncertainties. SPARE- 
ICE by Holl et al. (2014) showed a larger error in general, 
but their algorithm could quantitatively evaluate the contri- 
bution from each source, including CloudSat itself. The two 


algorithms are not directly comparable since we are under 
different metrics. Our algorithm meets the goal presented at 
the beginning, that to make cross-platform, cross-instrument 
consistent retrievals. 

5 Conclusions 

A fast empirical forward model built upon 7' cir -lWP relation- 
ships at MHS 157, 183.3 ± 3 and 190.3 GHz channels is de- 
veloped and used to retrieve tropical cloud IWP from MHS 
radiance measurements. The Tcir -IWP relationships at these 
channels are dependent on cloud top height h x in the tropics 
(Fig. 4), and the algorithm for retrieving IWP and lr t simulta- 
neously can improve the IWP accuracy. The IWP PDFs from 
MHS and CloudSat retrievals agree quite well, as expected 
for this constrained empirical forward model, over a wide dy- 
namic range of cloud ice (IWP = 0.5-10 kgm -2 , Fig. 6). The 
retrieval errors are also about the same magnitude (smaller 
than 100%, Fig. 10). The empirical forward model is valid 
for clouds with h t lower than 18 km and IWP greater than 
0.5 kgm -2 , but only in the tropics between 30° S and 30° N 
at present (Fig. 9). Beyond that latitude range, temperature 
lapse rate variations need to be taken into account to refine 
the 7’cir-lWP relationship. In addition, the algorithm is not 
accurate for retrieving IWP over elevated and arid topogra- 
phy (Fig. 8). 

Producing a CloudSat-consistent MHS IWP product has 
several direct benefits and important implications for study- 
ing clouds. Firstly, it helps to extend CloudSat cloud cov- 
erage with a wider swath width because frequent sampling 
from different operational satellites will allow frequent up- 
dates of fast-evolving weather phenomena such as hurri- 
canes and frontal systems. The new data can be used to im- 
prove weather prediction (e.g., cloud diurnal cycle) and long- 
term regional climate monitoring (e.g., IWP trend). Sec- 
ondly, our improved IWP retrieval method renders gener- 
ally larger IWP values than the NOAA operational product 
(Figs. 6 and 7). The approach we implemented with high- 
frequency microwave channels improves cloud detection in 
scenes with high IWP. Compared with CloudSat monthly cli- 
matology as well as the single -orbit measurements, we found 
that our results are closer to the CloudSat integrated ice wa- 
ter path. Thirdly, we show that replacing the 89 GHz chan- 
nel with 183.3 GHz channels for cloud ice retrieval reduces 
false detection of ice clouds and improves sensitivity to IWP 
as the higher-frequency channels are more sensitive to ice 
particle scattering. Lastly, the derived empirical T c i r — IWP 
relationships can be used to evaluate RTM simulations of 
cloudy-sky radiances, validate model assumptions, and im- 
prove model skills for data assimilation applications in the 
future (Fig. Cl). 

Although the empirical T c i r — IWP relationship developed 
here was from NOAA- 18 MHS, it is applicable to the sim- 
ilar channels used by other AMSU-B/MHS instruments on 
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the NOAA and Metop operational satellites for obtaining a 
longer data record and more frequent coverage. It can also 
be applied to other instruments that have the same combina- 
tion of channels, for example the Advanced Technology Mi- 
crowave Sounder (ATMS) onboard the Suomi-NPP satellite, 
or the Special Sensor Microwave Imager/So under (SSMI/S) 
onboard the Air Force F-16, F-17 and F-18 satellites. Fur- 


thermore, the approach we demonstrated in this study can be 
applied to IR/VIS sensors with the measurements collocated 
with CloudSat, such as Aqua Atmospheric Infrared Sounder 
(AIRS) and MODIS, to extend the sensitivity to lower IWP 
values and enhance the dynamic range of remote sensing of 
cloud ice from space. 
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Appendix A: Converting off-nadir Trir to nadir Tcir 

The collocation occurrences near outmost views are too few 
to form a statistically meaningful ^cir -IWP relationship, as 
shown in Fig. 2. An alternative way is to seek a conversion 
factor such that the PDF of converted off-nadir 7’ cn - matches 
that of the nadir Tcir, in which case the nadir 7’ cir -IWP rela- 
tionship can then be applied. 

If the cloud layer can be assumed to be plane -parallel and 
is not completely opaque to MHS, PDF s id e • cos f = PDF na dir 
roughly holds (thick solid lines in Fig. Ala), where t) is the 
zenith angle. This relationship means that, with equal proba- 
bility of seeing a cloud, the off-nadir cloud T c i r is colder than 
that from the nadir view, mainly due to the longer integration 
length of cloud water paths. However, such a relationship be- 
gins to lose its validity as the ice clouds become opaque at 
oblique views. As shown in Fig. Alb, the PDF curves be- 
come much more flat at T c „ = — 100 K compared with those 
at T C a = — 50 K. A factor of cos(£/2.1) can roughly capture 
such a view dependency. In the most extreme case, where 
7 C i r saturates at all view angles, the PDFs should become in- 
dependent of the view angle, where the factor can be set as 
cos(f /oo). 

By defining F = cos(£/fac), we can convert the off-nadir 
Tcir to the equivalent nadir Tcir by multiplying the former by 
F, i.e., 

Tcirlside ' F = T c ; r | equivalent nadir- (Al) 


0.0012 
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Figure Al. PDF of in the range of [—100, — 99 K] (a) and [—50, 
— 49K] (b) as a function of scan angle derived from a month of 
NOAA-18 MHS T c \ r data (December 2010). Thick solid curves are 
calculated from the mean PDF values averaged over the 10 nadir- 
view FOVs divided by a factor of cos(? /2.1) (a) and cost; (b), re- 
spectively, and are used to fit the observed PDF curves, f is the 
zenith angle. 




examination of other months of T cn PDFs for this channel 
shows the same warm bias (named as AT). At off-nadir 
views, AT becomes smaller, which follows the theoretical 
clear-sky limb-darkening curve AT na dir-cos£ = A7’ s ,a c . This 
means that the estimated Tccr for clear skies has a systematic 
error at 157 GHz, which could originate from the inaccuracy 
of MERRA surface emissivity or errors inside CRTM. To ac- 
count for this offset, all 7’ C]r values for CH#2 are subtracted 
from this offset before carrying out the retrieval: 


fac = f(T clr , CH#) is a function of T ca and the choice of 
channel number, which is estimated from the observation to 
follow such relationships: 


fac(CH#2) = 

' 0.7764 •e“ 0 00777 ’ dr 
1 

100 

fac(CH#4) = 

0.0013 ■ e“ 0 I0347dr 
• 1 
100 


if T c ir > -120 K 

if 0.7764 • e -°- 0077r cir < i 

if Tdr < -120 K 


if T cir > -120 K 
if 0.0013 • e _0 10347 i:ir < i 
if Tdr < -80 K 


(A2) 


(A3) 


fac(CH#5) = 

0.8160 •e“ 0 00987 ’ dr 
1 

100 


if Tcir > -120 K 


if 0.8160 • e -° 00987 cir < i (A4) 


if T cir — — 1 20 K 


Appendix B: Correction of T C i r for CH#2 

As shown by the black line in Fig. la, the PDF peak for Tcir 
of CH#2 is at 2 K instead of 0 K at the nadir view. Further 


Tcir(CH#2) — T c ir(CH#2) C alculated ~ 2 • COS£ (Bl) 


Appendix C: Comparison with RTM simulations 

As one of the most important motivations of this work, it is 
worthwhile checking whether the existing operational RTMs 
can reproduce the observed 7’cir-IWP relationship. Some re- 
search RTMs that have more sophisticated radiative trans- 
fer treatment produced comparable results with observations 
(e.g., Davis et al., 2007; Kulie et al., 2010). However, some 
of the major operational RTMs still have large biases in high- 
frequency microwave channels, which result in poor usage of 
cloudy/precipitating scenes observed by instruments such as 
MHS and SSMI. The goals of this exercise are hence to quan- 
tify operational model uncertainties and to determine how 
reliable these RTMs are for calculating cloudy-sky radiances 
as observed by satellite microwave channels/instruments. We 
would also like to identify possible causes of these model er- 
rors and develop a quick remedy for these operational mod- 
els. Since this problem by itself is very complicated, it will 
only be touched briefly in the appendices here, and a paper 
to explore further details is underway. 

Two popular models, CRTM and CRM, will serve as 
the representatives among different RTM variants. CRTM is 
widely used as the centerpiece of data assimilation in United 
States weather forecast systems, and CRM is used every day 
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to generate MLS cloud products. In both models, three ice 
clouds are fed in one by one, with cloud bottoms at 7.5 km 
and cloud tops at 10, 12 and 14 km, respectively. Two cloud 
shapes (convection with anvil cloud top; Gaussian shape) are 
tested. Since their results differ little, the first cloud shape 
is applied to all following studies. US standard atmosphere 
in the tropics is used as the background atmosphere. For 
the cloud droplet size distribution, CRM has several options. 
Only McFarquhar-Heymsfield (MH) and Gamma distribu- 
tions are tested with different combinations of parameter val- 
ues. MFl distribution was applied to deliver MLS IWC prod- 
ucts (Wu et al., 2009), and Gamma size distribution was as- 
sumed for CloudSat ice water product retrievals (Austin et 
al., 2009). In CRTM, only the cloud ice effective radius is 
tunable with fixed-width Gamma size distribution assump- 
tion. Both models and the observation are compared at nadir 
view only. 

Comparing 7’ cir responses from different channels to the 
same cloud is a straightforward yet very effective way of 
presenting many of the differences. As one can see from 
Fig. Cla, CH#2 should have a larger response to thin and 
medium thick clouds, and the penetration depths of CH#2 
and CH#5 are about the same when they encounter dense 
and thick clouds. Raw 7 b from the two channels showed the 
same features (not shown). CRTM produces almost identical 
reposes for CF1#2 and CH#5, both of which are however too 
weak compared with the observation (triangles). CRM pro- 
duces comparable dynamical ranges of 7' cir , with an effective 
radius of 160 pm and a width parameter of 2 (refer to Evans 
et al. (1998) for the format of the Gamma distribution), but it 
always generates a weaker response for CH#2, contradicting 
the observation (crosses), while the ?cir -IWP relationship for 
CH#5 is simulated quite well. The main caveat of CRM is 
that it only considers scattering while ignoring the ice emis- 
sion. The observed CH#2 / CH#5 T cn ratio suggests that ice 
cloud emission offsets as much as 30% of the cloud scat- 
tering impact for thin and medium thick clouds, while cloud 
scattering dominates for dense and thick clouds. Moreover, 
liquid droplets in mixed-phase clouds could contribute more 
to the emissions and further reduce the CH#5 T c ■„ response. 
That may explain the difference between CRTM and the ob- 
servation. 

The models simulate the observed 7cir ratio between CH#4 
and CH#5 better, as shown in Fig. Clb. Nevertheless, mod- 
els tend to over-predict the CH#4 response. Because CH#4 
is closer to the 183.3 GFlz water vapor absorption line, it is 


(a) <•» 




Figure Cl. Scatter plots of Tcir relationships between CH#2 and 
CH#5 (a), and CH#4 and CH#5 (b) from observed T c i r at MHS 
nadir view (black dots), simulated T c i r from CRTM (colored trian- 
gles) and from CRM (colored stars). Blue/cyan/red colors represent 
cloud layer with tops at 10, 12 and 14 km and bottoms at 7.5 km in 
the simulations. 


more sensitive to water vapor variations than CH#5. There- 
fore, the air above the cloud top might have been drier than 
the ambient air to lead to a smaller magnitude of CH#4 Tcir- 
This derivation is supported by some observational evidences 
for deep convective clouds (e.g., Chae et al., 2011). Again, 
CRTM overall produces weaker 7' c ,,-. probably due to the 
heavy extrapolation from thin cloud LUT to thicker clouds, 
and the narrow fixed-width Gamma size distribution. 

To summarize, the two widely used operational RTMs are 
far behind the empirical model in capturing the observed 
7’cir— IWP relationships for tropical ice clouds. For relatively 
thick ice clouds that should produce 30 K or more 7 b depres- 
sions, the warm bias of CRTM easily exceeds 100 %, making 
it impossible to assimilate the observed thick/precipitation 
ice clouds from high-frequency microwave instruments. For 
CRM, the biases can be as large as 50% for the 150 GHz 
channel. Some plausible explanations are given to explain 
the observed model discrepancies, which include but are not 
limited to the wrong extrapolation of the look-up table of ice 
clouds, oversight of emission from liquid droplets and ice 
particles, indifference to humidity above clouds and ambient 
air, and too narrow a width of Gamma size distribution. One 
should be advised that Gamma distribution is not indicated as 
the best cloud ice particle size distribution. Rather, it reflects 
the inter-model consistency between the two RTMs and the 
original assumptions made for CloudSat IWP retrievals. Due 
to the popularity of these two RTMs, this part of work is of 
much interest to a broad community. 
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